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Abstract 

A concise method for following the evolving geometry of a moving 
surface using Lagrangian coordinates is described. All computations can 
be done in the fixed geometry of the initial surface despite the evolving 
complexity of the moving surface. The method is applied to three prob- 
lems in nonlinear elasticity: the bulging of a thin plate under pressure 
(the original motivation for Foppl-von Karman theory), the buckling of a 
spherical shell under pressure, and the phenomenon of capillary wrinkles 
induced by surface tension in a thin film. In this last problem the inclu- 
sion of a gravitational potential energy term in the total energy improves 
the agreement with experiment. 



1 Introduction 

The elasticity theory of thin shells is largely differential geometry by another 
name. In this paper I describe a method for following the differential geometric 
data of a surface as it moves, and illustrate its application to non-linear elasticity 
theory. The equations of the method are completely general for smooth surfaces, 
and so could in principle describe the complex motions of crumpling up to the 
formation of singularities. 

Problems involving elastic membranes have been approached in several ways, 
including numerical simulation by triangulated surfaces, using a polyhedral ap- 
proximation to differential geometry pQ. Another approach has been to use 
differential geometry and scaling laws to understand the line and point singu- 
larities of crumpled surfaces analytically [2j |H [5j [6l [7] , and numerically [8]. 
The method of this paper generalizes familiar methods of mechanical engineer- 
ing for the non-linear elasticity theory of thin shells [HI [TU1 [H] in going beyond 
second order, and in treating initially curved surfaces in a unified way. 

Section [2] establishes notation for the differential geometry of a moving sur- 
face and shows how to use Lagrangian coordinates to simplify its description, 
the main idea of this paper. Section [3] summarizes the observations of the previ- 
ous section in a system of differential equations for the evolving surface and its 
strains. Section |4] compares this approach to Foppl-von Karman (FvK) theory, 
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and solves the motivating problem for that theory, the bulging of a thin rectan- 
gular plate subject to pressure, by integrating the evolution equations forward 
in time. Section [5] uses second order expansions of the crumpling equations to 
describe the buckling of a sphere under pressure. Section [5] uses the insights of 
Cerda and Mahadevan [T^] to give a more detailed description of a phenomenon 
recently discussed in [13j , capillary wrinkles induced by surface tension in a thin 
film. A previously unnoticed discrepancy with experiment is partially resolved 
with the inclusion of the gravitational potential energy of the system. 

2 Geometrical Methods 

In terms of smooth coordinates (ar, x 2 , x 3 ) in space one can describe the de- 
formation of a material object by the trajectories of its constituent particles, 
solutions of equations of motion 

where 

V = V*8i (2) 

is the vector field generating the flow, and t is a parameter along the flow. Inte- 
grating the system forward to t = 1, one can also think of V 1 as a displacement, 
a slight abuse of notation that should be clear from context. Metric relations 
among the particles are given by 

ds 2 = gijdx l dx-' (3) 

where </y is a Riemannian metric tensor, perhaps, but not necessarily, the Eu- 
clidean metric. 

I coordinatize the material object by Lagrangian coordinates, convected by 
the flow, i.e., every material point keeps the same coordinates that it had orig- 
inally. In this case the changing metric relationship of material points, namely 
the change in the expression Eq. ([3]) , is due entirely to the change in the metric 
components g^, because dx % , which for this purpose simply assigns to a line 
segment the coordinate difference of its endpoints, is invariant. The rate of 
change as a consequence of this deformation in the components of the metric g, 
or of any second rank tensor G, expressed in convected coordinates, is given by 
the Lie derivative|14j|15j 

£ v G(dj,d k ) = VG jk + G([d 3 ,V],d k ) + G(dj, [d k ,V]) (4) 

Here [ , ] is the Lie bracket of vector fields. It is more common to express 
objects like this, derivatives of tensors which are themselves tensors, in terms of 
the covariant derivative with respect to the metric connection, and to employ 
the conventions of raising and lowering indices with gij and its matrix inverse 
g tJ , such that, for example the covector with components 

Vt = 9ijVi (5) 
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has covariant derivative with respect to x k (denoted V$ ; fc), in terms of the ordi- 
nary partial derivative (denoted Vi jk ) given by 

V i .,k = V i , k +Ti k V j (6) 
where the coefficients of connection T are 

^ik = 2 gjm ^ %k ' m ~ gmi ' k ~ 9km,i) (7) 

It is straightforward to verify for any second rank tensor G^ v that 

£ v G(d k ,dt) = g l i {V 3 G m + Vj. ik G u + V^Gm) (8) 

In particular, if G is the metric tensor <?, which is a covariant constant, we 
recover the well known result 

£ v g(d k , di) = V t , k + V k , t = 2U M (9) 

where U is the rate of strain tensor of the flow V (or the first order strain of 
the displacement V). Nothing said above was specific to three dimensions, and 
therefore every statement can be interpreted as referring to a surface with a 
Riemannian structure if the indices take only two values and not three. From 
now on I shall use Latin indices for tensors in three-space, and Greek indices 
for tensors on a surface. 

Now consider a smooth material surface M, so thin that one may regard it 
as 2-dimensional, and let (x 2 , x 3 ) be coordinates in this surface, while x 1 = z is 
displacement along the normal to the surface, with the positive direction chosen 
conventionally, such that the surface M is z = 0. Such a coordinate system 
exists for a neighborhood of M such that \z\ < 1/C, where C is the supremum 
over M of both principal curvatures in absolute value. The metric tensor in 
these coordinates takes the form 

= ( 1 ° 

The tensor <? M „, with Greek indices taking values (2,3), is the first fundamental 
form of M, h^ v is the second fundamental form, and k^ v — hrjixv is the third 
fundamental form. All these tensors are associated with the surface M, and 
not with the ambient space. They do not depend on z, i.e. all z dependence in 
Eq. (|10p is explicit. The plus sign on the middle term is a conventional choice. 
On a sphere, for example, one could take the positive direction for z to be the 
outer normal direction, and the principal curvatures of the sphere to be positive. 

Now let a vector field (a, V^) be prescribed on M with normal component 
a (a; 2 , a: 3 ) and tangential components V^(x 2 ,x 3 ), and extend it to a neighbor- 
hood of M as 

W = ad z + V^dp - zG^a.^dv , (11) 

where initially the tensor G ^ = g^ u . In a short time At, the flow generated by 
the velocity field W changes the metric tensor components by approximately 

Ag = At£ w g (12) 
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The tensor g + Ag regarded as a tensor on 3-space expresses the ambient Eu- 
clidean geometry in Lagrangian coordinates. If g+ Ag is restricted to the surface 
z = and indices (2,3), one has the slightly altered first fundamental form of 
M 

G^v = (fiV + &9ii.v)\z=a - ( 13 ) 

expressing the non-Euclidean geometry of the slightly altered M induced by its 
embedding in the ambient Euclidean space. The term linear in z in Eq. (jlip 
was chosen to maintain the block diagonal form of Eq. (|10[) to first order in z. 
Therefore, taking the z-derivative, one has the slightly altered second funda- 
mental form of M, 

(14) 

The third fundamental form could not be computed in this way, but it is deter- 
mined by H^, 

= H^H\ . (15) 

I now imagine taking a sequence of such small steps, and I will continue to 
denote by G M „ and H^ v the evolving first and second fundamental forms giving 
the Riemannian structure on M induced by the embedding in Euclidean space. 
I will not make use of this Riemannian structure for computations, however. 

There is another natural Riemannian structure on M , namely that given by 
the original, undeformed first fundamental form (j^y-, together with its associated 
connection, etc., which I shall continue to use, being careful not to give it erro- 
neous interpretations. This Riemannian structure, unlike G^ u , has no obvious 
geometrical meaning on the deformed surface, but it is still useful in a formal 
way. Another possible interpretation, deliberately suppressing the geometrical 
meaning of G M „ , is to imagine a surface that is not deformed by the flow W but 
carries tensor fields G^ v and H^ v , initially coinciding with g M „ and h^ v , that 
are deformed by W . That these happen to be the first and second fundamental 
forms of an evolving surface is forgotten. In this picture the undeformed g^ v has 
an obvious geometrical meaning as the metric on the underlying, unchanging 
surface which is the arena for the evolving G Mt , and H^. 

In Eq. (TTTj) I introduced the tensor G A " y , initially g^ u . More generally G ^ 
is the inverse of G^ as a matrix. It is a tensor field on M, but it is not obtained 
from G^y by raising indices. Raising indices is an operation accomplished by 
g^, my chosen Riemannian structure, not by C". The prime on G' is a 
reminder that it is not some version of the tensor G. 

I have shown how G^ and change, to first order, under a deformation 
(a, V M ) of M, assumed now always to be extended off M as in Eq (jlip . In 
turn, (a, V^) might evolve so as to reduce at each step a free energy functional 
depending on G^ and H^. In this way I will arrive at crumpling equations, 
a system of differential equations for (a, V 7 *), G^, and H^, describing the 
evolution of M. Before considering the equation for (a, V^), though, there is 
another issue to consider. 



-(g^ + Ag^) 
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This formulation leaves implicit what the evolving surface actually looks 
like, since mere knowledge of and is not a convenient description of 
M. To keep track of the positions of points on the surface, one should integrate 
Eq. fl]) using components of W(0, x^) — (a, V^) with respect to fixed Cartesian 
coordinate axes. Let X A (x 1 , x 2 , x 3 , t) be a Cartesian coordinate function in 
space. It is time independent in the physical sense, but its functional form 
depends on time because the x l evolve in time. The 1-form dX A = X A ^dx^ 
assigns the X A component W A to the vector W . This 1-form evolves in time 
at the rate given by the Lie derivative 

£ w dX A (di) = WdX A (di) + dX A ([d l ,W}) (16) 
= W s X A tij + U Y \ (17) 
= (X^WSli (18) 

Thus U A = U l X A , the X^-component of any vector field U = U l di at time t, 
can be found using X A (x x , x 2 ,x 3 , t) solving 

with appropriate initial conditions. By the definition of the coordinate x 1 = z, 
the Cartesian coordinate X A is an affine linear function of z. It is essential 
therefore to expand X A W^ only to first order in z in Eq. fTS}. To be completely 
explicit, X A is independent of z and we can represent 

X A = Y A {x\x*)+zZ A {x\x*). (20) 



Then Eq. US} says 



dXf 

dt 

dY A 
dt 

dZ A 
dt 



= Z A V^ - Y A a^G " u (21) 
= (X A a + X A Vn^ (22) 
= {Z A V V - Y A a,xG'" x )^ (23) 



The linear approximation I have made in the neighborhood of M obscures the 
fact that if W were made to carry affine normal lines to affine normal lines 
exactly, as one could always require by a suitable nonlinear extension W of 
(a, V^) off M, then X A W^ would be exactly an affine linear function of z 
without approximation. The evolution of M is the same for any extension, 
however, so what looks like a linear approximation in the method is actually 
exact. 

As a special case, I describe motion at constant velocity, i.e., dW A /dt = 
for each component A. Then 

-A, 



^jxfp^ mwl + xf s^ (24) 
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Thus the components of W must evolve according to 



= -XliX+WtliW* (25) 

Here X k A is the inverse of Xj 4 , considered as a matrix. Eq. (f2l)|) for straight line 
motion is recognizable as 

<9W fc 

+ W 3 ViW k = (26) 

dt 

where is the covariant derivative with respect to the metric connection of the 
Euclidean metric in 3-space expressed in the evolving Lagrangian coordinates. 
I emphasize that I have chosen, however, not to use the evolving geometry but 
rather the fixed initial geometry of M for all computations, a great simplifica- 
tion. 



3 Evolution Equations 

By the arguments of the previous section the surface M evolves according to 
dC 

= WG KXili + V%G ltX + V%G lut +2aH KX (27) 

OH \ 1 / 

— a K K \ — a a ; k + 2°"'^ ^" (~ G K \- V + G U \ K + G VK -\) 

+VH nX . li + V%H„ X + V^H^ (28) 

(29) 

Using these relations one can find how other geometric quantities change, for 
example the area element \fG dx 2 dx 3 , involving the determinant of the first 
fundamental form 

G = G22G33 — G23G32 (30) 

The result is 

^ = (V»VG),» + aG'^H^VG (31) 
Integrating one finds \fG and hence dilation strain. The strain tensor 

\(G^-g^) (32) 



can be found by integrating Eq. (|27jl . A natural definition for nonlinear shear 
strain S^ LV is 



dS^ = 1 dG^_ 1 dVG 
dt ~ 2 { dt VG dt 



G^ . (33) 



The subtracted term removes the contribution of dilation strain. is not 



traceless, in general, beyond first order. 



G 



4 Comparison with Foppl-von Karman Approach 



A simple example illustrates the use of this formalism and points out its rela- 
tionship to Foppl-von Karman (FvK) theory [9] . FvK considers the equilibrium 
state of a thin membrane subject to external forces and boundary conditions. 
Since the metric strain within a membrane is typically small, even for large 
normal displacements, it makes sense to continue to use linear stress-strain re- 
lationships. The strain may, however, be a nonlinear function of displacement, 
and therefore displacement may be nonlinearly related to stress. FvK thus 
produces nonlinear equations for the equilibrium shape of an elastic membrane 
subject to external stress. 

Historically this idea was implemented by expanding the strain tensor to 
first order in tangential displacement but second order in normal displacement. 
I derive the FvK strain by solving the evolution equations to first order in 
and second order in a, continuing to use the notation of previous sections, with 
the initial velocity vector 

= ad z + V^^df, - zG^a^d v (34) 

of Eq. (fTTj) . I am using the superscript W to indicate the initial value, which is 
also the zeroth approximation for an iterative solution. Other initial values are 
9nu = G$ = S^ v and = H/S) = 0. I use Picard's method to generate the 
solution to the differential system Eqs. (TIT))) , (|125l ([27]) , and (|28[) iteratively as a 
power series in t, taking M to be the Euclidean plane with the usual Cartesian 
coordinates. In this case there is no distinction between indices up and indices 
down, and covariant derivatives are ordinary partial derivatives. Iterating once, 
and ignoring quadratic terms except in a gives 

= S KX +t(V^+V^ K ) (35) 

= -**M (36) 
Vf> = taa !tl (37) 

Iterating a second time, still ignoring quadratic terms except in a, gives 

Finally, evaluating at t = 1, gives the FvK metric strain 

\{G { 2 - S KX ) = \{V^ X + + a<K a, x ) - a kA (39) 

This is the computational starting point for FvK theory. The rest of that theory 
follows from minimizing the elastic energy, expressed as a quadratic functional 
of this strain and the first order bending strain H^tJ , to find the equilibrium 
shape. 

The approach of this paper is to develop the nonlinear strain as the solution 
to a differential system. From that point of view the derivation of Eq. |39|) is not 
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very natural, since to obtain it one must artificially impose the condition that 
the trajectories of the particles are straight lines, a condition that introduces, 
via Eq. (|3~T|) , a second order correction into the strain that is necessary to obtain 
Eq. (f3U)) . Although one can certainly parameterize the possible final shapes of 
M by displacement of particles along straight lines, it is a different thing to say 
that particles actually move along straight lines. FvK theory does not claim 
this, and in that sense it is not a dynamical theory. A dynamical theory would 
determine the evolution of the velocity vector (a, V^) by some local physical 
law, replacing Eq. ([25]) in the differential system. It would be a simpler theory, 
both conceptually and computationally, in that solving it would only require 
integrating a differential system forward in time. I will do the obvious thing 
and choose W to reduce the elastic energy at each step, seeking the minimum. 
A typical phenomenological elastic energy functional is 

E = E d + E s + E c (40) 

where 

e - - K,(f <4i) 

E s = n [ S KX S KX ^gdx 2 dx 3 (42) 
Jm 

E c = ~J (d^H^-gThp) 2 y/gda?dx* (43) 

and where A, ji, k are the 2D compression modulus, shear modulus, and bending 
modulus respectively. The area element involves y/g, not \/G, because the 
energy due to metric strain is better understood to be per unit mass, not per 
unit area, and the mass is convected with the material coordinates. The system 
will move, if possible, to lower its energy, so one must compute the variation of 
E with respect to a small normal displacement Sa and tangential displacement 

~ 5E „ „ SE 



5E= / — SV + -Z- Sa J-gdx z dx A (44) 
Jm YoV^ da J 

The work done on M in deforming it represents energy given up by some other 

part of the system, so this work should be added with a minus sign to the total 

change in energy. Work done by pressure P in a small normal deformation Sa, 

for instance, is 

W = p[ SaVGdx 2 dx 3 (45) 

where now one must use the physical area element VGdx 2 dx 3 on M. A small 
displacement in the direction opposite to this "gradient," i.e. 

V - (47) 
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will lower the energy and move the system toward a local minimum. The linear 
operators L a and Ly include a projection onto the space of admissible veloc- 
ity vector fields. They must define positive semi-definite quadratic forms with 
respect to the inner product given by integration over M. Apart from these re- 
quirements, they will vary with the application. This is just the familiar notion 
of conjugate gradient. One could also think of L a and Ly together as defin- 
ing a generalized mobility tensor, because it transforms generalized force into 
velocity If one only wants to know the final state, one could try to choose L a 
and Ly so as to reach equilibrium in the most efficient way. In any case, the 
dynamics of the system is not completely determined by the elastic energies, 
and additional physical considerations must be added to complete the theory in 
a specific application. 

Eqs. (|46|) and (|47|) . together with the evolution equations of Section [3l are 
what I mean by Lagrangian crumpling equations. The original problem ad- 
dressed by FvK theory, the bulging of a square plate fixed on the boundary 
and subject to pressure, can be solved straightforwardly in this way. Represent 
all geometric data by discretization on a square grid of points of the original 
square. Spectral methods (fast Fourier transform with anti-aliasing) make the 
computation efficient, and the gradient flow converges quickly to a solution. 

5 Buckling of a sphere under pressure 

I consider an elastic spherical shell subject to pressure P, described by the 
phenomenological energies of Eqs. (|¥T|) . P2")) . (|4*3")) . and (|4"5|) . For small enough 
pressure the sphere is uniformly compressed, but as pressure increases it buckles. 
I will describe the buckling by using expansions of strain to second order in 
displacement where necessary, not the FvK expansion, but the "dynamic" one 
of this paper, found by solving the crumpling equations iteratively. It turns out 
that the expansion must include more terms than FvK. 

For a sphere of radius R, in terms of spherical polar coordinates (#,</>), 



Taking R = 1, and regarding all quantities now as dimensionless, the perturbed 
geometric quantities in a general displacement (a, V*) are 




diag(i? 2 ,i? 2 sin 2 6) 
9^/R 



(48) 
(49) 
(50) 



H, 



G, 



g/iu + Vn; V + V V ;n + lag^v 

+ + V vtli 



(51) 
(52) 



Vg 



y/g[l + (V* + 2a) + -{V%V% + WV^ 
AaV^ + 2V»a ilL - aa^ + 2a 2 )] 



+ 



(53) 



The area element \fG had to be found to second order in displacement. To first 
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order in displacement the shear strain in the sphere is 

S^^liV^ + V^-V^g^) (54) 
Parameterize the displacement by coefficients (af m ,i< m ,C| m ), such that 

a = ^2a em Y tm (55) 

lm lm 

where the Yg m are spherical harmonics and 632 = —£23 = sin#, £22 = £33 = is 
the antisymmetric tensor. Then for example the change in the mean curvature 
of the perturbed sphere is 

SH = G^H^ - gTh^ = + - 2 K™^m (57) 

lm 

so that the curvature energy is 

E * = \ E[^+ 1 )-2] 2 |^™| 2 (58) 

lm 

It vanishes for I = 1, as it must by Galilean invariance, and it is independent of 
the tangential displacement V^. The other energy expressions are 

W = 4TTPa oY o-PY,£(e+l)at m bi m + 2PY,\a em \ 2 

ira £m 

£(e+l)\a em \ 2 + 2\a lm \ 2 } (59) 

lm 



Ed = |^H(£+l)6£ m + 2o £m ] 2 



2 

lm 



+ Aa 0Q Y 0Q J2l-^(£ + l)aimbim + i(l + l)\aim\ 2 + 2\a £m \ 2 ] (60) 



lm 



E s = |^^+l)[^+l)-2](|6 £m | 2 + | Qm | 2 ) (61) 

lm 

These expansions have been carried out to second order in all coefficients, but 
they anticipate that a no is the same order as |a^ m | 2 for I > 1, so that some 
terms quadratic in aoo appear to be third order. I also anticipate that the first 
response to pressure is a uniform compression 

|ooo| ~ j (62) 
so that consistency requires P << A. Now seek the minimum of the total energy 

Etot = W + E d + E s + E c (63) 
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by choice of (ae m , be m , Q™). Ignoring corrections of order P/A gives 

*m = (64) 

2Aft£ m 

lm = M{l+l) + n[l(i+\)-2] (65) 



A 

£m 



1 



4/x 



A^+l) + ^(f+l)-2] 



(66) 



Putting these expressions back into £tot gives 

£tot = -ttP 2 /A + (67) 

where Q is a diagonal quadratic form in the coefficients ai m . One must deter- 
mine the sign of the diagonal elements in Q, since the appearance of negative 
coefficients in Q corresponds to the onset of buckling in the corresponding mode 
I. Introducing the notation 

a = e(£+l), (68) 
the diagonal element is (a — 2)F(a), where 

, N P 2A y u 2 (a-2) + 2^oA 2 , 

= "2 + [A. + ^-2)p + < a - ^ <™ 

It is clear that for any a > 2 the diagonal element becomes negative for large 
enough pressure P, so that buckling must occur, but the only relevant value of 
a is the one for which this first happens as P increases. If k > (k) ct , where 

(«)« = ^ (l + f ) (70) 

(I have restored dimensional factors of R), then F(a) is monotonically increasing 
for a > 2. Thus as P increases, F(a) first becomes negative for the lowest 
nontrivial shape mode I = 2 corresponding to a = 6, and the buckling will be 
of quadrupole shape. If, on the other hand, k < (k) ci , the more interesting 
case, then F has a local minimum for some a > 2, and hence a buckling mode 
that doesn't simply grow from the translation mode but appears at a higher a. 
Values of (a, P) for which F(a) has a double root correspond to the onset of 
buckling into this mode. Solving F(a) — and F'(a) = simultaneously, and 
restoring dimensional factors R, one finds the buckling mode lb and buckling 
pressure Pb 

2 M « + y/2i?2 /iK A(A + fi) 



ni?3 = V32^(A + M )-4^ m4bW ^ (?2) 
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According to Eq. (|TTj) , the wavelength A of the buckling mode has the form 
A ~ (k/A) 1 / 4 argued by Cerda and Mahadevan [12], although the mechanism is 
not quite the same as the one they describe. In their case the applied stress is 
anisotropic, while here the symmetry breaking is spontaneous. 

The second order expansion does not determine how the crumpling pro- 
ceeds once buckling has occurred, but it does give an initial condition for the 
crumpling equations, which are now just a differential system of equations for 
(aim,b£ m ,C£ m ). Solving this system numerically might be tractable, since the 
the right hand side of the system involves only integrals of products of spherical 
harmonics and their covariant derivatives over the sphere. 

In doing this problem of the buckling sphere, I noticed that second order 
expansions of strains on curved surfaces must include terms not only quadratic 
in normal displacement a, which is the FvK prescription, but also mixed terms 
like aV^. Eq. (I53|) contains such a term, for example, since nothing is omitted 
there, through second order. If such terms are mistakenly ignored, second order 
expansions of elastic energies fail to be Galilean invariant for the simple reason 
that in translating a curved surface normal and tangential displacements are 
necessarily of the same order. It is not true that tangential displacements are 
small even when normal displacements are large, which is the FvK argument 
for ignoring them. This problem with the translation mode (£ = 1) also affects 
nearby fs, by continuity. That the expansion is then only accurate for large I 
means that it is good only if the wavelength of the perturbation Ye m is much 
less than the radius of curvature, but this is just the case in which we can regard 
the surface as flat. That is, simply generalizing the FvK strain of Eq. (|39|) to a 
curved surface is not much of an advance over assuming the surface to be flat. 

6 Capillary Wrinkles 

A recent paper described radial wrinkles produced in a floating thin film by the 
surface tension of a small drop placed at the center [13] . Through a combination 
of physical arguments, dimensional analysis, and experiment, the phenomenon 
was explained in a sufficiently quantitative way to become a useful assay for the 
properties of the film. The theory is that of Cerda and Mahadevan [12]. In that 
paper, physical intuition simplifies the problem, which is essentially a problem 
of FvK theory, but at the expense of some of the details. The methods of this 
paper, guided by the intuition of [T5], stay closer to FvK theory and show in a 
little more detail how the result emerges. I also incorporate a term that turns 
out to be important but that was not included in the original discussion, the 
gravitational potential energy of the supporting fluid, disturbed by the wrinkling 
film. 

A thin film disk of radius R floats on a water surface, subject to surface 
tension a, and a small water drop of radius p is placed at its center. It is 
equivalent to think of an annular film subject to radial stress a at its outer radius 
R and radial stress 2a at its inner radius p. Choose units so that p = 1, and use 
cylindrical polar coordinates {z, r, 9) in the same formalism as in other sections. 
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The equilibrium state is attained by a displacement (a, V^) that minimizes the 
total energy E, given by the sum of the elastic energies, Eqs. (|4Tj) . (|42j) . and 
(T33|) , the negative of the work done by surface tension 



/-2tT /-2tT 

W = 2a V 2 { Pl 9)pd9-a V 2 (R,9)Rd9, (73) 
Jo Jo 



r 2ir r-2-K 

V 2 (p,9)pd9 - a \ 
'o Jo 

and the gravitational potential energy of the water 



2 r2ir nR 

E g = - J J Pw ga 2 rdrd9 (74) 

where pw is the density of water. 

Expand the total energy E in the sense of FvK theory, that is, use strains 
linear in V* and quadratic in a. Then in the absence of wrinkling (i.e., a = Q), 
E is minimized by the radial displacement 

V 2 = -- + Br (75) 
r 

with 

A ~ 2pjR? - 1) ' B ~ 2A(i?2 _ i) ( 7b > 

in which the disk is slightly dilated and sheared. In the process the energy is 
lowered by 

2 M (i? 2 -l) 2A(i? 2 -l) ' 1 ' 

A still lower energy is attained, however (and this is a variational estimate), by 
a state with m radial wrinkles of the form 

acos(m6>) 

a - —p (78) 

o A m 2 a 2 , , 

v - -~ + Br -WTW^ (79) 

3 ma 2 sin(2m6>) 

V ^2P ( 80 ) 

The form of V^ 1 is chosen to cancel the m 2 term in the shear strain due to the 
wrinkle a (that is the idea of [12] translated into the language of this paper). 
One is still free to choose the parameters A, B, (3, to, and a, this last being 
the dimensionless amplitude of the wrinkles. The dependence on A and B is 
quadratic, so that the best values are trivially found. The energy E then has 
the form 

2 

E=^2 E i3 (l3)a 2i m 2j (81) 

i,J=0 
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with Eqi = Eq2 = 0. Minimizing with respect to a 2 leads to 



2 E n + 2E 12 m 2 



a ~ E 21 +2E 22 m2 m 

Since the coefficient E\ 2 in the numerator comes from the bending energy alone, 
it is negligible compared to En, and thus in the wrinkling regime a ~ 1/m, in 
agreement with the intuition of [12] . Substituting this value back into Eq. (|8Tj) , 
one finds that the optimal m 2 satisfies a cubic equation 

= A im 6 + A 2 m 4 + A 3 m 2 + A 4 (83) 

Finally one should seek the optimal value for (3. For all physically reasonable 
values of the parameters in the problem, the optimal value turns out to be (3 = 
(only approachable as a limit) corresponding to wrinkles that keep a constant 
amplitude. 

The optimal values computed above turn out to be insensitive to the dilation 
strain, due to equal stress a at inner and outer radii, and only sensitive to the 
unbalanced stress a in the center due to the drop, creating shear strain, the 
only strain that can be relieved by wrinkling. Taking the limit p — > 0, and also 
ignoring k/A and n/fi, since the bending modulus is small for thin films, leads 
to simple values for the coefficients in Eq. 



Ax = -2E 12 E 22 « -tt 2 Are/64 (84) 

A 2 = -3E 12 E 21 « -3ttVk/256 (85) 

A 3 = 2E 22 E W - 4E 20 E 12 - E U E 21 w tt 2 p w gR 2 A/64 + n 2 [icj/256 (86) 

A 4 = E w E 21 -2E 20 E 11 ^7r 2 Pw gR 2 t i/256 (87) 

The roots of Eq. ([55]) . for typical physical values, are determined almost entirely 
by A\ and A3, that is, 

m^i-A./A,) 1 / 4 (88) 

in formal agreement with |12) . There are two regimes, depending on the rel- 
ative importance of the gravitational term in A3. If the gravitational term is 
unimportant, 

/ u \l/4 /C r\i/4 

771 £H I I 



V4A/ 



Since \x ~ A, the dimensionless first factor is about \j\/2 « 0.7. This factor was 
measured experimentally in [13] and found to be about 3.6. Because of the 4th 
root, the discrepancy is very large. If the gravitational term dominates in A3, 



m 



f^ R2 \ lf \^ (90 ) 



with a crossover between the two regimes at 



R ~ \ l T— ~ 1-4 mm, (91) 
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taking the value a = 72 x 10 3 J/m 2 from [T3]. Since R in that experiment 
was 11.4 mm, it was in the second regime, and the dimensionless first factor in 
Eq. is roughly 2, still not 3.6, but closer! 

Like reference [12] , this analysis does not explaine the observed length of the 
wrinkles, which seem here to have length R. A solution going beyond second 
order, obtained by integrating the crumpling equations forward in time, might 
resolve this question. 
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